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1 Introduction 



During the last decade, a quantum tunneling of magnetization has become 
among the most actively studied topics in a condensed matter physics. The 
immense interest aimed at better understanding of this quantum phenomenon 
has been mainly stimulated by a recent experimental observation of the quan- 
tum spin tunneling in a large number of single-molecule magnets (see Ref. 
[1] and references therein). By the term single- molecule magnet one denotes 
an assembly of weakly interacting clusters of magnetic metal atoms that usu- 
ally possess an extraordinary strong magnetic anisotropy. Hence, the single- 
molecule magnets often provide very good examples of magnetic systems with 
a strong uniaxial magnetocrystalline anisotropy, i.e. so-called Ising-like spin 
systems. Of course, the Ising-type interaction by itself cannot be a source of 
the quantum spin tunneling experimentally observed in these systems. It turns 
out, however, that this quantum phenomenon arises in the most cases due to 
the higher-order crystal-field terms. According to a number of experimental 
and theoretical studies it is now quite well established that the spin tunneling 
observed in Fe 4 [2], Fe 8 [3], Fei 9 [4], or Mn 4 [5] clusters originates to a major 
extent from a second-order biaxial crystal-field potential. 

Extensive studies focused on the magnetic properties of small clusters shed 
light on the effect of single-ion anisotropy terms D (uniaxial anisotropy) and 
E (biaxial, also called rhombic anisotropy). In contrast to a quite well un- 
derstood role of both single-ion anisotropies in the small magnetic clusters 
(zero-dimensional systems), the situation becomes much more complex and 
also obscure in one- and two-dimensional spin systems. In fact, ground-state 
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properties of a spin-S* Ising model with the rhombic crystal-field potential E 
have been only recently examined by Oitmaa and von Brasch within an effec- 
tive mapping to the transverse Ising model [6]. On the basis of this effective 
mapping, a zero-temperature quantum critical point can be exactly located 
in the one-dimensional model, while for the two-dimensional models, they can 
be estimated with a high numerical accuracy using the linked-cluster expan- 
sion method [6,7]. Nevertheless, the magnetic behavior of these models has 
not been investigated at non-zero temperatures beyond the standard mean- 
field and effective-field theories [8], random phase approximation [9], or above 
mentioned linked cluster expansion [10]. It should be stressed that the biax- 
ial anisotropy essentially influences magnetic properties of a large number of 
polymeric molecular-based magnetic materials, too. From the most obvious ex- 
amples one could mention: NiF 2 [11], NiN0 3 .6H 2 [12], Ni(CH 3 COO) 2 .4H 2 
[13], Mn(CH 3 COO) 2 .3H 2 [14], CoF 2 [15], CoCl 2 .6H 2 [16] and a series of 
compounds Fe(dc) 2 X [17], where X stands for halides and dc for the dithio- 
carbamate or diselenocarbamate groups, respectively. 

In this article, we will focus on the effect of uniaxial and biaxial crystal-field 
potentials affecting the magnetic behavior of mixed spin- 1/2 and spin-3/2 hon- 
eycomb lattice. When assuming the Ising-type exchange interaction between 
nearest neighbors, the model becomes exactly solvable within an extended 
star-triangle mapping transformation. Thus, the considered model provides 
a noble example of statistical system, which enables to study an interplay 
between quantum effects and temperature in a spontaneously ordered mag- 
netic system. Moreover, a magnetic structure of the mixed-spin honeycomb 
lattice occurs rather frequently also in the molecular magnetism, what clearly 
demonstrates a large family of polymeric two-dimensional compounds of fol- 
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lowing chemical formula: A I M n M III (C 2 04)3 [18], where A 1 stands for a non- 
magnetic univalent cation N(C n H 2n +i)4 or P(C n H 2n +i)4 (n = 3 — 5), M 11 and 
M m denote two- and three-valent metal atoms Cu 11 ^ = 1/2), Ni 11 ^ = 1), 
Co 11 ^ = 3/2), Fe n (,S = 2) or Mn n (,S = 5/2) and respectively, Cr 111 ^ = 3/2) 
or Fe III (5' = 5/2). Actually, it turns out that the crystal structure of these 
polymeric molecular-based magnetic materials consists of the well-separated 
two-dimensional layers, in which regularly alternating M n and M m magnetic 
metal atoms constitute more or less regular honeycomb lattice (Fig. 1). In 
consequence of a large magnetocrystalline anisotropy of these materials, one 
should also expect a relatively strong uniaxial (Ising-like) anisotropy, as it has 
already been suggested in the theoretical studies based on the effective-field 
theory and Monte-Carlo simulations [19]. Hence, the magnetic compounds 
from the family of oxalates A I M II M III (C20 4 )3 represent good candidates to be 
described by the proposed model. 

The outline of this paper is as follows. In the next section, a detailed descrip- 
tion of the model system is presented and then, some basic aspects of the 
transformation method will be shown. Section 3 deals with a physical inter- 
pretation of the most interesting results and finally, some concluding remarks 
are drawn in Section 4. 

2 Model and method 

Let us consider a magnetic structure of the mixed-spin honeycomb lattice, as it 
is schematically depicted in Fig. 1. To ensure exact tractability of the model 
system, we will further suppose that the sites of sublattice A are occupied 
by the spin-1/2 atoms (depicted as full circles), in contrast to the sites of 
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sublattice B that are occupied by the spin-3/2 atoms (open circles). Assuming 
the Ising-type exchange interaction J between nearest-neighboring spin pairs 
only the total Hamiltonian of the system reads: 

3N N N 

n = jJ2 sm + dj: {sif + ej: i(s£f - (sin (i) 

(k,j) k£B k&B 

where N is a total number of sites at each sublattice, £1* and S%(a = x,y,z) 
denote standard spatial components of the spin-1/2 and spin-3/2 operators, re- 
spectively. The first summation in Eq. (1) is carried out over nearest-neighboring 
spin pairs, while the other two summations run over the sites of sublattice B. 
Apparently, the last two terms D and E are the crystal-field potentials that 
measure a strength of uniaxial and biaxial anisotropies acting on the spin- 
3/2 atoms. It is also worthy to note that there is one-to-one correspondence 
between the Hamiltonian (1) and the effective spin Hamiltonian with three 
different single-ion anisotropy terms D x , D y and D z : 

3N N N N 

n = JT, $m + DZ H (^) 2 + DX H (st) 2 + Dy £ (sf) 2 . (2) 

(k,j) k&B k&B k&B 

As a matter of fact, one can easily prove the equivalence between (1) and (2) 
using following mapping relations between the relevant interaction parameters: 

D = D z — -(D x + D y ), and E = -{D x - D v ). (3) 

It should be also mentioned here that by neglecting the biaxial anisotropy, 
i.e. setting E = in Eq. (1) or equivalently D x = D y in Eq. (2), our model 
reduces to the exactly soluble model settled by Gongalves [20] several years 
ago. Accordingly, the main attention will be focused here on the effect of bi- 
axial anisotropy, which influences thermodynamical and dynamical properties 
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in a crucial manner. Really, the £'-term related to the biaxial crystal-field 
anisotropy should cause non-trivial quantum effects, since it introduces the x 
and y components of spin operators into the Hamiltonian (1). As a result, it 
is responsible for the onset of local quantum fluctuations that are obviously 
missing in the Ising model with the uniaxial crystal-field potential D only. 

It is therefore of interest to discuss an origin of the biaxial crystal-field po- 
tential E. The origin of this anisotropy term consists in the low-symmetry 
crystal field of ligands from the local neighborhood of the spin-3/2 atoms. 
It is noticeable that a threefold symmetry axis oriented perpendicular to the 
honeycomb layer prevents an appearance of biaxial crystal-field potential in a 
regular honeycomb lattice with a perfect arrangement of the oxalato groups, 
as well as magnetic metal atoms. However, a small lattice distortion which oc- 
curs rather frequently in the low-dimensional polymeric compounds due to the 
Jahn- Teller effect can potentially lower the local symmetry. In consequence of 
that, the distortion of lattice parameters can be regarded as a possible source 
of the biaxial crystal-field anisotropy. The most obvious example, where the 
lattice distortion removes the threefold symmetry axis represents the single- 
molecule magnet Fe4, in which three outer Fe atoms occupy two non-equivalent 
positions around one central Fe atom [2] . 

Let us turn our attention to the main points of the transformation method, 
which enables an exact treatment of the model system. Firstly, it is very 
convenient to write the total Hamiltonian (1) as a sum of site Hamiltonians: 

N 

n = J2 n (k) , (4) 
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where the each site Hamiltonian fi^ involves all interaction terms associated 
with one spin-3/2 atom residing on the kth site of sublattice B: 



WW = S' h E k + {S z k fD + [(^) 2 - (Sl) 2 }E, 



(5) 



with E k = J(fi kl + jl%2 + friz)- While the Hamiltonians (5) at different sites 
commute with respect to each other {\f&\iiS^\ = 0, for each i ^ j), the 
partition function of the system can be partially factorized and consequently 
rewritten in the form: 



N 



k=l 



(*)1 



(6) 



In above, (3 = 1/(/cbT), k B is Boltzmann's constant, T the absolute temper- 
ature, Tr^} means a trace over spin degrees of freedom of sublattice A and 
Trs k stands for a trace over spin states of the kth spin from sublattice B. 
A crucial step in our procedure represents the calculation of the expression 
Tr 5fc exp[— (3H^]. With regard to this, let us write the site Hamiltonian (5) 
in an usual matrix representation: 



V 



9D 3E k 



o 

V3E 




\ 



D E k 
— + — 
4 2 





V3E 



V3E 


D _Ek 
4 2 





V3E 




9D 3E k 
U 4 2 J 



(7) 



in a standard basis of functions | ±3/2), | ±1/2) corresponding, respectively, to 
the four possible spin states S k = ±3/2, ±1/2 of the kth atom from sublattice 
B. Although it is easy to find eigenvalues of the site Hamiltonian (7), with 
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respect to further calculation it is more favorable to obtain directly matrix 
elements of the expression exp[— fiiiS k \ When adopting the Cauchy integral 
formula, one readily attains the matrix elements for an arbitrary exponential 
function of the site Hamiltonian (7): 



A 



. .1 



(exp{aH ik) ] 

A u = exp[a(jD + ±E k )] (cosh«+) + sinh«+)}, 



A 2 2 = exp 
A 33 = exp 
A 44 = exp 



4" 



-sinh(a4)}, 



a(-D --E k )] {coshK,-) 



\D + l -E k )] {cosh«+) - sinh«+)}, 



a 1 



D--E k )]{cosh(a^) + 
5 „ 1 „ ,i y/3E . 



D-E k 



sinh« fe )}, 



A 13 = A 31 = exp[a(-L> + -E k ) 
A 2 4 = A42 = exp [a(-D - -E k ) 



sinh(a^), 
VSE . 



— sinh« fc ), 



(8) 



where £ k = ^ (D ± E k ) 2 + 3E 2 and a marks an arbitrary function. After 
substituting a = —(3 in the set of Eqs. (8), the relevant trace Tr 5fc exp[— (3H^] 
can easily be calculated. Moreover, its explicit form immediately implies a 
possibility of performing the standard star-triangle mapping transformation: 



Tr Sfe exp[-/3W (fc) ] = 2exp[-5/3 J D/4 - (3E k /2\ cosh^p + E k ) 2 + 3E 2 
+ 2 exp[-5/3£>/4 + [3E k /2] cosk(j3y/{D - E k ) 2 + 3E 2 ) = 



= A exp[(3R(fi z kl fx z k2 + /4 2 /4 3 + fx z k3 fx z kl ) 



(9) 



which replaces the partition function of a star (i.e. the four-spin cluster con- 
sisting of one central spin-3/2 atom and its three nearest-neighboring spin-1/2 
atoms) by the partition function of a triangle (i.e. the three-spin cluster com- 
prising of three outer spin-1/2 atoms in the corners of equilateral triangle), 
as shown in Fig. 1. A physical meaning of the mapping (9) is to remove all 
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interaction parameters associated with the central spin-3/2 atom and to re- 
place them by the effective interaction R between the outer spin-1/2 atoms. It 
is noteworthy that both mapping parameters A and R are "self-consistently" 
given by the transformation equation (9), which must be valid for any com- 
bination of spin states of three spin-1/2 atoms. In consequence of that, the 
transformation parameters A and R can be expressed as: 

A = /Jfl = h(|), do) 

where the functions $1 and $ 2 are defined as follows: 



$i=2 exp(-5/3JJ>/4 - 3/3 J/4) cosh(/ty(3 J/2 + D) 2 + 3E 2 
+ 2 exp(-5/3£>/4 + 3(3 J/4) cosh(/3y/(3J/2 - D) 2 + 3E 2 ) , 

$2 = 2 exp(-5/3/J/4 - (3 J/4) cosh J/2 + D) 2 + 3E 2 ) 

+ 2 exp(-5/3/J/4 + ^ J/4) cosh^^ J/2 - /J) 2 + 3E 2 ) . (11) 

When the mapping (9) is performed at each site of the sublattice B, the 
original mixed-spin honeycomb lattice is mapped onto the spin-1/2 Ising tri- 
angular lattice with the effective interaction R given by the "self-consistency" 
condition (10)-(11). Indeed, the substitution of equation (9) into the partition 
function (6) establishes the relationship: 

Z(P,J,D,E) = A N Z t (P,R), (12) 

between the partition function Z of the mixed-spin honeycomb lattice and the 
partition function Z t of the corresponding spin-1/2 triangular lattice. Above 
equation constitutes the basic result of our calculation, since it enables rela- 
tively simple derivation of all required quantities, such as the magnetization, 
quadrupolar moment, correlation function, internal energy, specific heat, etc. 
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In addition, by combining (12) with (9) one readily proves a validity of fol- 
lowing exact spin identities: 



(fM, /4> •••> /4)> = (fM, v z j, -,»k))t, (is) 

(f2(S k , S\, S z k , /j, z kl , n z k2 , fj, z k3 )) = 

I TrsMSjj, SI S z k , fjfojjfa fa) expj-^W] 
\ Tr 5fc expf-^W] 

where (...) represents the standard canonical average over the ensemble defined 
by the Hamiltonian (1) and (...) t the one performed on the spin-1/2 Ising 
triangular lattice with the effective exchange interaction R. Here, f\ is an 
arbitrary function of the spin variables belonging to the sublattice A, while 
f 2 denotes an arbitrary function depending on the kth spin of sublattice B 
and its three nearest neighbors from the sublattice A. By applying the spin 
identity (13), one straightforwardly attains the following results: 



m A = (Afci) = (£fci)t = m t , (15) 

CA = (AfelAfc2> = (AfclAL)* = <h, (16) 
tA = (AfciALAD = (frtl^P-Wt = U, (17) 

whereas the second spin identity (14) enables after some algebra derivation of 
quantities depending on the spin variable from the sublattice B: 

m B = (S z k ) = -3m A (#i + K 2 )/2 - 2t A {K l - 3K 2 ), (18) 
77 ee ((S z k ) 2 ) = (K 3 + 3K 4 )/4 + 3c A (K 3 - K A ). (19) 

In above, ttia (%) labels the single-site magnetization of the sublattice A 



(B), 7] denotes the quadrupolar moment and finally, Ca and t A stand, respec- 
tively, for the static pair and triplet correlation functions between the relevant 
spins of sublattice A. Obviously, the exact solution of both sublattice magne- 
tization and quadrupolar moment require only a knowledge of the single-site 
magnetization mt, nearest-neighbour pair correlation function q and triplet 
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correlation function t t on the corresponding spin-1/2 Ising triangular lattice 
unambiguously given by R (10)-(11). Fortunately, the appropriate exact so- 
lution of these quantities is well known, hence, one can directly utilize final 
results derived in Refs. [21]. Finally, the coefficients emerging in the previous 
set of Eqs. (17)- (20) are listed below: 

K 1 = F 1 (3 J/2 + D, 3 J/2 - D, 3 J/2), K 2 = F 1 (J/2 + D, J/2 - D, J/2), 
K 3 = F 2 (3J/2 + D, 3 J/2 - D, 3 J/2), K 4 = F 2 (J/2 + D, J/2 - D, J/2), 

(20) 

where the functions Fi(x,y,z) and F 2 (x,y,z) are defined as follows: 



x sinh(/Va: 2 + 3E 2 ) 

Fi(x,y,z)- 



+ 



F 2 ( x,y,z) 



Vx 2 + 3E 2 cosh(f3Vx 2 + 3E 2 ) + exp((3z) cosh(/V?/ 2 + 3E 2 ) 
y sinh(/Vy 2 + 3E 2 ) 

vV + 3E 2 exp(-/3z) cosh(/Vz 2 + 3E 2 ) + cosh^^ 2 + 3E 2 ) 

1 cosh(/Vz 2 + 3E 2 ) - exp(pz) cosh(/Vl/ 2 + 3E 2 ) _ 

2 cosh(/5\/a?TW) + exp((3z) cosh(/3y/y 2 + 3E 2 ) ' 
y sinh(/3Vy 2 + 3E 2 ) 



vV + 3E 2 exp(-/3z) cosh(/Vx 2 + 3E 2 ) + cosh^vV + 3E 2 ) 



^5 x sinh(/3 v / x 2 + 3£ 2 ) 



4 Va; 2 + 3£ 2 cosh(/5v / x 2 TW) + exp(/?z) cosh(/5 v ^ J TW) ' 

(21) 



Now, we will derive an exact result for one dynamical quantity, namely, the 
time-dependent autocorrelation function. It should be noted here that exactly 
tractable models offer only seldom a possibility to investigate their spin dy- 
namics. On the other hand, the dynamical quantities such as autocorrelation 
and correlation functions are important also from the experimental point of 
view, because their magnitude directly determines a scattering cross section 
measured in inelastic neutron scattering experiments [22], or a spin- lattice 
relaxation rate provided by a nuclear magnetic resonance technique [23]. 
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As a starting point for calculation of the time-dependent autocorrelation func- 
tion Cau to can, for convenience, serve the exact spin identity (14): 



1 /Tr Sh {[S' k (0)Sm + S* k (t)S' k (0)\ exp[-/W( fc )]}\ 



2\ ' T.vcxp W«] """ /' 
where the symmetrized form in the definition of C** to is used to construct 
a Hermitian operator, S k (t) = exp(^ iL )S'| exp(— ^ iL ) represents the Heisen- 
berg picture for the time-dependent operator S k (t), h stands for the reduced 
Planck's constant and % = \J — 1. Next, the matrix elements of expressions 
exp(±^ iL ) can be in turn evaluated by putting a = ±^ into the set of Eqs. 
(8). Then, after a straightforward but a little bit tedious calculation, one ar- 
rives to a final result for the dynamical autocorrelation function: 



C :i t0 = \[K 5 (D + 3 J/2, D - 3J/2, 3 J/2) + 3K 6 (D + J/2, D - J/2, J/2)] 
+ 3ct[K 5 (D + 3 J/2, D - 3J/2, 3 J/2) - K 6 (D + J/2, D - J/2, J/2)], (23) 



where the time-dependent coefficients K 5 and K 6 are, for brevity, explicitly 
given in the Appendix. 



3 Results and discussion 



Before proceeding to a discussion of the most interesting results, it is notice- 
able that the results derived in the previous section are rather general, since 
they are valid for a ferromagnetic (J < 0) as well as ferrimagnetic (J > 0) 
version of the model under investigation. In what follows, we shall restrict 
our analysis to the ferrimagnetic model only, since the polymeric compounds 
from the family of oxalates [18] fall mostly into the class of ferrimagnets. Nev- 
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ertheless, it appears worthwhile to remark that a magnetic behavior of the 
ferrimagnetic system completely resembles the one of ferromagnetic system. 
Finally, one should also emphasize that the mapping (9) remains invariant un- 
der the transformation E <-> — E. For this reason, one may consider without 
loss of generality the parameter E > and consequently, x-, y- and z-axis 
then represent the hard-, medium- and easy-axis for a given system. 

3. 1 Ground-state properties 

At first, we will take a closer look at the ground-state behaviour. Taking into 
account the zero-temperature limit (T — > + ), one finds the following condition 
for a first-order phase transition line separating two different magnetically 
ordered phases denoted as OPi and OP 2 : 



Moreover, one easily attains from Eqs. (15)-(21) analytical results for the 
single-site sublattice magnetization (m A , m B ), the total magnetization nor- 
malized per one magnetic atom m = (rriA + m B )/2 and the quadrupolar 
moment r], as well: 




(24) 



OPi: m A 



1 




3 
2 



V 
J 



(25) 




m 



1 



3 
2 



V 
J 



7) 



5 
4 



+ 



3 
2 



D 
J 





OP 2 : m A 



1 




(26) 
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For better illustration, Fig. 2 depicts the ground-state phase diagram in the 
E-D plane (Fig. 2a) and the zero-temperature variations of the magnetization 
and quadrupolar moment with the biaxial anisotropy when D/J = 1.0 (Fig. 
2b). It should be mentioned that by neglecting the biaxial anisotropy, i.e. 
by putting E/J = 0.0 into the phase boundary condition (24), one recovers 
the boundary uniaxial anisotropy D/J = 0.75 in accordance with the results 
reported by Goncalves [20] several years ago. Moreover, the OPi (OP2) phase 
corresponds in this limit to the simple ferrimagnetic (antiferromagnetic) phase 
with both sublattice magnetization oriented antiparallel with respect to each 
other: m\ = —0.5, m# = 1.5 in the OPi and respectively, uia = —0.5, m# = 
0.5 in the OP 2 . Apparently, the spin-3/2 atoms occupy exclusively the | + 3/2) 
(1 + 1/2)) state in the OPi (OP 2 ) phase when E = is satisfied. 

The situation becomes much more complex by turning on the biaxial anisotropy. 
Even though the sublattice magnetization uia remains in the ground-state at 
its saturation value in both OPi and OP 2 phases, the sublattice magnetization 
rriB is gradually suppressed by the effect of biaxial anisotropy (see Fig. 2b). It 
is quite obvious from this figure that the biaxial anisotropy gradually destroys 
the perfect ferrimagnetic (antiferromagnetic) spin arrangement, which occurs 
in the OPi (OP 2 ) phase in the limit of vanishing E. Let us find a primary 
occasion for this unexpected behavior accompanied by a spin reduction at 
sublattice B. According to Eq. (26), one finds rj — m B = 3/4 to be valid in 
the whole parameter space corresponding to the OPi. From an elementary 
consideration it can be easily understood that the spin-3/2 atoms must oc- 
cupy in the OPi phase either the | + 3/2) or | — 1/2) state in order to satisfy 
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simultaneously both values of and r] in the ground state. Aforementioned 
argument is also supported by the fact that the E-teim. does not couple in 
the Hamiltonian (7) the | +3/2) and | — 1/2) states with the | + 1/2) and 
| — 3/2) ones. Hence, the | + 3/2) and | — 1/2) states are in the OPi phase 
the only allowable states, while the | + 1/2) and | — 3/2) states are, on the 
contrary, the forbidden ones. The observed spin reduction at sublattice B can 
be thus attributed to the local quantum fluctuations induced by the biaxial 
anisotropy, which in turn lead to a spin tunneling between the | +3/2) and 
| — 1/2) states in the OPi. Although the occupation of the minority | — 1/2) 
state rises steadily with increasing the biaxial anisotropy strength, it is no- 
ticeable that the | +3/2) state still remains the state with the most probable 
occupation. For completeness, it should be pointed out that in accord with our 
expectation the negative uniaxial anisotropy (D < 0) reduces the occupation 
of the minority | — 1/2) state, while the positive one (D > 0) favors it. 

Quite similar situation emerges also in the OP2 phase. However, it is worth- 
while to remark that the OP2 phase appears in the region of strong uniaxial 
anisotropies D / J > 0.75 only. Due to a strong positive uniaxial anisotropy, 
the spin-3/2 atoms undergo a well-known spin transition from the | +3/2) to 
I + 1/2) state, which macroscopically manifests itself in the phase transition 
from the OPi to OP 2 phase. As stated before, the biaxial anisotropy couples 
together the | + 1/2) and | — 3/2) states and therefore, the tunneling between 
these spin states should be expected to occur in the OP 2 . The analytical so- 
lution for rriB and rj (27), as well as a validity of the relation rj + rriB = 3/4 in 
a whole parameter space corresponding to the OP 2 phase, indeed confirm this 
suggestion. However, the negative (positive) uniaxial anisotropy prefers (re- 
duces) the occupation of minority | — 3/2) state in the OP 2 phase in contrast 
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to the abovementioned trends observed in the OPi phase. 

Now, let us step forward to the discussion of the time dependent autocorrela- 
tion function (23). Among other matters, this quantity can serve in evidence 
whether the spin-3/2 atoms fluctuate in the OPi (OP2) phase between their 
allowable | + 3/2) and | — 1/2) (| + 1/2) and | — 3/2)) spin states. Unfor- 
tunately, it is quite tedious to derive from Eq. (23) a simple analytical ex- 
pression for C^uto i n the zero-temperature limit, hence, we report for C** to 
numerical results obtained at very low temperature (k B T / J = 0.001) close to 
the ground state. Fig. 3 displays the time variation of C** to for several values 
of uniaxial and biaxial crystal-field potentials. Since C%* t0 evidently varies in 
time, it clearly demonstrates the zero-temperature spin dynamics between the 
allowable states. Moreover, a detailed analysis reveals that C** to is in the zero- 
temperature limit a harmonic function of time and whence, the time depen- 
dence can be characterized by an angular frequency u± = ^-\J (-j ± 1) + 3(-j) 2 
depending on whether the system resides the OPi, or OP2 phase. This result is 
taken to mean that the spin system necessarily recovers after some character- 
istic recurrence time (r± = 2it/uj±) its initial state. More specifically, Fig. 3a 
illustrates the time variation of C** t0 in the OPi phase, because with respect 
to Eq. (24) one never approaches the OP 2 phase for D/ J = 0.0. Fig. 3a clearly 
clarifies the role of biaxial anisotropy: the stronger the biaxial anisotropy E/J, 
the greater the angular frequency of spin tunneling and in the consequence of 
that, the shorter the appropriate recurrence time. Furthermore, the increasing 
strength of the biaxial anisotropy enhances also the amplitude of oscillation in 
the time-dependence of C%* to . This observation would suggest that the increase 
of biaxial anisotropy enhances a number of the spin-3/2 atoms tunneling dur- 
ing the recurrence time between the | + 3/2) and | — 1/2) states in the OPi 
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phase. However, since the equilibrium magnetization does not change in time, 
the number of atoms that tunnel from | +3/2) to | — 1/2) state must definitely 
be the same as those that tunnel from the | — 1/2) to | + 3/2) state. 

To illustrate the effect of uniaxial anisotropy, the time variation of C** t0 is 
shown in Fig. 3b for E/J = 0.5 and several values of D/J. Apparently, the 
Cauto oscillates for strong negative (positive) uniaxial constants D/J in the 
vicinity of boundary values C** to = 2.25(0.25). These values clearly demon- 
strate that the prevailing number of spin-3/2 atoms occupy in the OPi (OP 2 ) 
phase the | + 3/2) (| + 1/2)) state, since C** to = i] when t — 0. Moreover, 
the stronger the uniaxial anisotropy (independently of its sign), the smaller 
the relevant amplitudes of oscillation, i.e. the smaller the number of tunneling 
atoms during the recurrence time. On the other hand, the increasing strength 
of uniaxial anisotropy enhances the angular frequency of oscillation, what 
means, that the tunneling atoms return from the minority | — 1/2) (| — 3/2)) 
state to the most probable occupied | +3/2) (| + 1/2)) state of the OPi (OP2) 
phase after a shorter recurrence time. 

3.2 Finite-temperature behaviour 

In this part, we would like to comment on the finite-temperature behaviour 
of the system under investigation. Let us begin by considering the effect of 
uniaxial and biaxial single-ion anisotropies on the critical behaviour. For this 
purpose, two typical finite-temperature phase diagrams are illustrated in Fig. 
4a and 4b. In both figures, the OPi (OP 2 ) phase can be located below the phase 
boundaries depicted as solid (dashed) lines, while above these boundary lines 
the usual paramagnetic phase becomes stable. Further, open circles represent 



17 



special critical points at which both OPi and OP 2 phases coexist. Actually, 
we have not found any phase transition between the OPi and OP2 phases at 
non-zero temperatures, what indicates that the OPi phase coexists with the 
OP 2 one at non-zero temperatures merely for the same D / 'J — E / J values as 
in the ground state (24). Finally, a closer mathematical analysis reveals that 
both temperature-driven phase transitions which are related to the OPi and 
OP 2 phase, respectively, are of a second order and belong to a standard Ising 
universality class. 

Fig. 4a shows the critical temperature (T c ) as a function of the uniaxial 
anisotropy for several values of the biaxial anisotropy. The critical temper- 
ature versus uniaxial anisotropy dependence is quite obvious for E/J = 0.0, 
when increasing D/J, T c monotonically decreases. Hence, the critical temper- 
ature approaches in the limit D/J — > — oo(+oo) its minimum (maximum) 
value ksTc/ J = 0.3796... (1.1389...) in agreement with the exact T c (triple T c ) 
of the spin- 1/2 Ising honeycomb lattice. A gradual decline of the transition 
temperature can obviously be explained as a consequence of the fact, that the 
positive uniaxial anisotropy energetically favors the low-spin | ± 1/2) states 
before the high-spin | ±3/2) ones. The most interesting finding to emerge 
here is that the biaxial anisotropy may significantly modify the critical be- 
havior of the studied system. As a matter of fact, T c firstly reaches its local 
minimum at certain positive D/J and then rises steadily to its limiting value 
k B T c / J = 0.3796.... The extraordinary increase of T c in the region D/J > 1.0 
can be explained through a suppression of the occupation of minority | — 3/2) 
state, which appears in the OP 2 phase due to the uniaxial anisotropy effect. In 
accordance with previous assumption, the greater the biaxial anisotropy (i.e. 
the greater a number of atoms that occupy the minority | — 3/2) state), the 
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more impressive increase of T c can be observed. In addition, it is easy to un- 
derstand from here that the biaxial anisotropy substantially lowers the critical 
temperature of OPi phase in the D/J < 0.0 region, in that it is responsible 
for the quantum spin tunneling between the | +3/2) and | — 1/2) states. 

To illustrate the influence of biaxial anisotropy on the criticality, the depen- 
dence of transition temperature on the biaxial anisotropy is shown in Fig. 
4b for several values of the uniaxial anisotropy D/J. As one would expect, T c 
gradually decreases with increasing the biaxial anisotropy for any D/J < 0.75. 
It is quite obvious that the appropriate suppression of T c can be attributed 
to the quantum fluctuations, which become the stronger, the greater the ratio 
E/J. Apart from this rather trivial finding, one also observes here the peculiar 
dependences with the non-monotonic behavior of T c . The critical temperature 
may exhibit only a slight variation with increasing E/J (as it is in the case 
of D/J = 1.0), or it may show unexpected local minima, as it is in the case 
of D/J = 2.0 and 3.0. Since the local minima can be located very near to 
the coexistence point of the OPi and OP 2 phases (depicted as open circles), 
the relevant increase of T c can be related to the OP2 — > OPi phase transition. 
Namely, the most populated | + 1/2) spin state in the OP2 is replaced after 
this phase transition by the | + 3/2) state, which is the most occupied spin 
state in the OPi. The spin crossover from the low-spin | + 1/2) to the high-spin 
I + 3/2) state must lead, of course, to a slight increase of T c . 

At this stage, let us provide an independent check of the critical behavior by 
studying thermal dependences of magnetization. The single-site magnetization 
is plotted against temperature in Fig. 5 for the biaxial anisotropy E/J = 0.5 
and several values of the uniaxial anisotropy D/J. Fig. 5a shows a typical sit- 
uation observed in the OPi phase: the more positive the uniaxial crystal-field 
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potential D/J, the stronger the spin reduction (the lower the magnetization 
TOb) due to the | + 3/2) <-> | — 1/2) spin tunneling. In consequence of that, 
the total magnetization alters from a standard Q-type dependence observed 
for D < (see for instance the curves for D/J = —2.0 and 0.0) to a more 
interesting R-type dependence, which occurs for positive uniaxial anisotropies 
(D/ J — 0.5 and 0.75). Unusual slope in the thermal dependence of total mag- 
netization can be related to a more rapid thermal variation of mg. In fact, on 
account of the quantum fluctuations m B is thermally easier disturbed than 
which, on the contrary, always exhibits the standard Q-type behavior (spin- 
1/2 atoms are not directly affected by the biaxial crystal-field potential E). 
Furthermore, Fig. 5b shows how the situation changes by considering the tran- 
sition toward the OP 2 phase. Actually, both magnetically ordered phases OPi 
and OP 2 have the same internal energy (coexist together) at D/J = y/l3/A 
when Ej J = 0.5, while the OP 2 phase becomes more stable if Dj J > \^13/A. 
Accordingly, Fig. 4b displays the thermal variation of sublattice magnetization 
exactly at the OP1-OP2 phase boundary and in the OP2 phase (D/J = 1.0 
and 1.5). The corresponding thermal dependences of total magnetization are 
plotted in the insert of Fig. 5b. As it is apparent from these figures, the initial 
value of rriB is suppressed from its saturation value (me = 0.5) owing to a 
presence of the minority | — 3/2) state. Nevertheless, a large number of spins 
can be thermally excited to the | +3/2) state for D/J from the vicinity of 
OP1-OP2 phase boundary and hence, rapidly increases upon heating (see 
the curve for D/J = 1.0). As a result of this thermal excitation, the total 
magnetization exhibits N-type dependence with one compensation point in 
which rriA and completely cancel out (see the insert in Fig. 5b). Finally, 
even for stronger uniaxial anisotropies (e.g. D/J = 1.5) the total magnetiza- 
tion recovers the Q-shape, since the thermal fluctuation prefer excitations to 
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the | — 1/2) state rather than to the | + 3/2) one. Such a thermal excitations 
must, naturally, lower the sublattice magnetization rag. 

To conclude our discussion devoted to finite-temperature properties, let us 
proceed to the time variation of dynamical autocorrelation function C** t0 as 
depicted in Fig. 6 for E/J = 0.5 and three selected values of D/J. In order 
to enable a comparison between the displayed data at various D/J, the rel- 
evant temperatures are normalized with respect to their appropriate critical 
temperatures, i.e. we have defined the dimensionless temperature r = T/T c 
that measures a difference from the critical point (r c = 1.0). The time varia- 
tions of C** t0 from Fig. 6a and 6b display the relevant changes of dynamical 
autocorrelation function in the OPi phase, while Fig. 6c shows the correspond- 
ing dependences in the OP 2 phase. It can be easily understood that C** to is 
not in general the time-periodic function at non-zero temperatures no mat- 
ter whether considering C** t0 in the OPi, or OP 2 phase. In fact, C** t0 arises 
according to Eq. (23) as a superposition of two harmonic oscillations with 
two different angular frequencies cu± = ^ \J~(jf± l) 2 + 3(-j) 2 and also various 
amplitudes. The interference between these harmonic oscillations gives rise to 
a rather complex time variation of C%* to , which is in general aperiodic, dis- 
playing nodes and other typical interference effects. The periodicity of C** to 
at non-zero temperatures is maintained only for some particular E/J — D/J 
values, which retain the ratio u + /lj_ to be rational, while in any other case, 
C*uto behaves aperiodically. 

The dependences drawn in Fig. 6 nicely illustrate also the temperature effect 
on the spin dynamics. It follows from these dependences that some ampli- 
tudes are suppressed as the temperature increases, while another ones become 
more robust. Obviously, in the high-temperature regime that amplitudes be- 
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come dominant, which coincide to the oscillation with lower angular frequency. 
Contrary to this, the amplitudes arising from higher frequency oscillation dom- 
inate in the low-temperature regime. The most miscellaneous time variation 
of C** t0 thus emerges in the vicinity of critical temperature (r 1.0), which 
represents an intermediate temperature range between the low- and high- 
temperature regime. However, a rather exceptional case is displayed in Fig. 
6c, where the most miscellaneous dependence appears surprisingly at sub- 
stantially lower temperature (r = 0.25) rather than the critical one (r c = 1.0). 
When looking back to the thermal variation of magnetization depicted in Fig. 
5b, one finds a feasible explanation for this striking behavior. It turns out 
that the temperature (r 0.25) of the most miscellaneous time variation of 
Cauto coincides with the temperature k B T '/ J ps 0.1, at which the most robust 
spin excitation to the | + 3/2) can be observed. In addition to the allowable 
| + l/2) and | — 3/2) states, a large number of the spin-3/2 atoms is therefore 
thermally excited to the | + 3/2) spin state. This observation would suggest 
that the thermal excitations can basically modify the spin dynamics as well. 

4 Concluding remarks 

In this article, the exact solution of the mixed spin-1/2 and spin-3/2 Ising 
model on honeycomb lattice is presented and discussed in detail. The particu- 
lar attention has been focused on the effect of uniaxial and biaxial crystal-field 
anisotropies acting on the spin-3/2 atoms. As it has been shown, a presence 
of the biaxial anisotropy significantly modifies the magnetic behavior of the 
system under investigation. It turns out that already a small amount of the 
biaxial anisotropy raises a non-trivial spin dynamics and basically influences 
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the thermodynamic properties, as well. 

The most striking finding to emerge here constitutes an exact evidence of 
the spin tunneling between the | ±3/2) and | =p 1/2) states in two different 
magnetically ordered phases OPi and OP 2 , respectively. Macroscopically, the 
tunneling effect decreases the critical temperature of the magnetically ordered 
phases and appreciably suppresses the magnetization of spin-3/2 sublattice. 
This quantum reduction appears apparently due to the local quantum fluctu- 
ations arising from the biaxial crystal-field potential. 

There is an interesting correspondence between the model described by the 
Hamiltonian (1) and a similar model with a local transverse magnetic field 
Q acting on the spin-3/2 atoms only [24]. However, similarity in their actual 
properties is not accidental, in fact, when neglecting the uniaxial crystal- field 
potential D in Hamiltonian (1), an effective mapping E <-> Vt ensures the 
equivalence between both the models. Since this mapping is not related to 
the magnetic structure in any fashion, the appropriate correspondence can be 
extended to the several lattice models. It is therefore valuable to mention that 
magnetic properties of the models with a local transverse field become a sub- 
ject matter of many theoretical works [25]. Apparently, the magnetic behavior 
of these systems should completely resemble that one of their counterparts 
with the biaxial crystal-field potential. 

Finally, let us turn back to the origin of biaxial anisotropy. Uprise of this 
anisotropy term in the mixed-spin honeycomb lattice is, namely, closely as- 
sociated with at least a small lattice distortion. To simplify the situation, 
the proposed Hamiltonian (1) accounts for the biaxial crystal- field anisotropy, 
while a difference between exchange interactions in the different spatial direc- 
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tions has been for simplicity omitted. Nevertheless, the developed procedure 
can be rather straightforwardly generalized to an anisotropic model account- 
ing for the different interactions along various spatial directions. Moreover, 
the biaxial anisotropy can be even considered as an arbitrary function (linear, 
quadratic, exponential, logarithmic, ...) of the ratio between appropriate in- 
teraction parameters. Hence, it would be very interesting to find out whether 
such a system is instable toward the spontaneous lattice distortion caused by 
the spin-Peierls phenomenon. In this direction continues our next work. 



5 Appendix 



An explicit form of the coefficients K 5 and K 6 is given by: 



K 5 (x,y,z) = G(x,y,z)/^ 1 , and K 6 (x, y, z) = G(x, y, z)/<b 2 , 



where the function G(x, y, z) is defined as follows: 



G(x, y, z) = ^[W^H^x, z) + W^H^y, -z)\ + ^[W 2 (x)H 2 (x, z) 
+ W 2 (y)H 2 (y, -z)\ - [W 3 (x)H 3 (x, z) + W 3 (y)H 3 (y, -z)\, 
W^x) = |x 2 + E 2 [l + 2 cos(2t v / a: 2 + 3E 2 /H)] }/(x 2 + 3E 2 ), 

W 2 (x) = {x 2 - 3E 2 [l - 2 cos(2tVx 2 + 3E 2 /H)] }/(x 2 + 3E 2 ), 



x 1 — cos (2tVx 2 + 3E 2 /h)] |/(x 2 + 3E 2 ), 
H x {x, y) = exp(-5/3D/4 - (3y/2) [cosh^Va; 2 + 3E 2 ) 



- x smh((3Vx 2 + 3E 2 ) /Vx 2 + 3E 2 \ , 
H 2 (x, y) = exp(-5/3L>/4 - (5y/2) [cosh^Vx 2 + 3E 2 ) 

+ x sinh(/V:r 2 + 3E 2 )/Vx 2 + 3E 2 ] , 
H 3 (x, y) = exp(-5(3D/4 - (3y/2)V3E sinh^Va; 2 + 3E 2 )/Vx 2 + 3E 2 . 
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Figure captions 

Fig. 1 The segment of a mixed-spin honeycomb lattice. The lattice positions of the 
spin- 1/2 (spin-3/2) atoms are schematically designated by the full (open) circles, 
the solid lines label interactions between nearest neighbors. The dashed lines 
represent the effective interactions between three outer spin-1/2 atoms arising 
after performing the mapping (9) at kth site. 

Fig. 2 a) The ground-state phase diagram in the E/J — D/J plane; b) The single-site 
sublattice magnetization \mA\ (dotted line) and tub (dashed line), the total single- 
site magnetization \m\ (solid line) and the quadrupolar momentum rj (solid line) 
as a function of the biaxial anisotropy E/J at T = and D/J = 1.0. 

Fig. 3 The time variation of dynamical autocorrelation function C^ to at very low tem- 
perature (ksT/J = 0.001) close to the ground state: a) for D/J = 0.0 and various 
E/J; b) for E/J = 0.5 and various D/J; Time axis is scaled in K/J units. 

Fig. 4 a) The dependence of critical temperature on the uniaxial anisotropy D/J for 
several values of biaxial anisotropies E/J; b) The dependence of critical temper- 
ature on the biaxial anisotropy E/J for several values of uniaxial anisotropies 
D/J. Solid (broken) lines correspond to critical temperatures of the OPi (OP2) 
phase. Open circles denote the critical temperatures for a such particular case, 
when both the ordered phases OPi and OP2 coexist in the ground state (see the 
text). 

Fig. 5 The thermal dependences of single-site magnetization for E/J = 0.5 and: a) 
D/J = -2.0, 0.0, 0.5 and 0.75; b) D/J = v 7 !^, 1.0 and 1.5. The dotted (dashed) 
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lines stand for the sublattice magnetization \mA\ ("Ib), the solid lines for the 
total single-site magnetization \m\. Fig. 5b shows the temperature dependences 
of sublattice magnetization |m^| and mg only, the insert shows the appropriate 
changes of total magnetization \m\. 

6 The time variation of dynamical autocorrelation function C** to when Ej J = 0.5 is 
fixed and D/J changes: D/J = —2.0 (upper), 0.0 (central) and 1.0 (lower panel). 
The relevant time variations of C^ to are displayed at three different temperatures, 
which are normalized with respect to their critical temperatures in order to get 
the ratio r = 0.75, 1.0 and 1.25 in Fig. 6ab and respectively, r = 0.1, 0.25 and 
1.0 in Fig. 6c. Time axis is scaled in the h/J unit. 
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Fig. 1 Strecka et al. 





Fig. 3 Jascur et al. 




Fig. 4 Jascur et al. 
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